#if !defined(IMGLIB_RASTER_H)
#define IMGLIB_RASTER_H

#include "typelib/types.h"

namespace std { struct bidirectional_iterator_tag; }

namespace img {

enum order_type { bottom_up, top_down };

inline order_type flip(order_type o)
	{ return order_type(1-o); }

typedef size_t size_type;

}

#include "imglib/operations.h"

#include <vector>
#include <functional>

namespace img {

template <class T, class E> class raster_iterator {
private:
	T *			_raster;
	size_type	_x, _y;
public:
	typedef E	value_type;
	typedef E 	pointer; //???
	typedef E 	reference; //???
	typedef ptrdiff_t	difference_type;
	typedef std::bidirectional_iterator_tag	iterator_category;
public:
	raster_iterator() {}
	raster_iterator(T & t, size_type x, size_type y)
		: _raster(&t), _x(x), _y(y) { }
	bool operator == (const raster_iterator & rhs) const
		{ return (_raster == rhs._raster && _x == rhs._x && _y == rhs._y); }
	bool operator != (const raster_iterator & rhs) const
		{ return (_raster != rhs._raster || _x != rhs._x || _y != rhs._y); }
	raster_iterator & operator ++() { 
		++_x; 
		if (_x == _raster->width()) 
			{ ++_y; _x = 0; } 
		return * this; 
	}
	raster_iterator & operator --() {
		if (_x == 0)
			{ --_y; _x = _raster->width(); }
		--_x;
		return * this;
	}
	E operator * () 
		{ return (*_raster)(_x, _y); }
	E rel(int dx, int dy) 
		{ return (*_raster)(_x+dx,_y+dy); }
public:
	T & raster()
		{ return *_raster; }
	const T & raster() const
		{ return *_raster; }
	size_type x() const
		{ return _x; }
	size_type y() const
		{ return _y; }
	operator types::point<size_type>() const
		{ return types::point<size_type>(_x, _y); }
	template <class U>
	types::point<U> point() const
		{ return types::point<U>(_x, _y); }
};


template <class T, class E, class RI = raster_iterator<T,E> > 
class const_raster_iterator {
private:
	const T *	_raster;
	size_type	_x, _y;
public:
	typedef E	value_type;
	typedef E 	pointer; //???
	typedef E 	reference; //???
	typedef ptrdiff_t	difference_type;
	typedef std::bidirectional_iterator_tag	iterator_category;
public:
	const_raster_iterator() {}
	const_raster_iterator(const T & t, size_type x, size_type y)
		: _raster(&t), _x(x), _y(y) { }
	const_raster_iterator(const RI & i)
		: _raster(&i.raster()), _x(i.x()), _y(i.y()) { }
	bool operator == (const const_raster_iterator & rhs) const
		{ return (_raster == rhs._raster && _x == rhs._x && _y == rhs._y); }
	bool operator != (const const_raster_iterator & rhs) const
		{ return (_raster != rhs._raster || _x != rhs._x || _y != rhs._y); }
	const_raster_iterator & operator ++() { 
		++_x; 
		if (_x == _raster->width()) 
			{ ++_y; _x = 0; } 
		return * this; 
	}			               
	const_raster_iterator & operator --() {
		if (_x == 0)
			{ --_y; _x = _raster->width(); }
		--_x;
		return * this;
	}
	E operator * () const 
		{ return (*_raster)(_x, _y); }
	E rel(int dx, int dy) const
		{ return (*_raster)(_x+dx,_y+dy); }
public:
	const T & raster() const
		{ return *_raster; }
	size_type x() const
		{ return _x; }
	size_type y() const
		{ return _y; }
	operator types::point<size_type>() const
		{ return types::point<size_type>(_x, _y); }
	template <class U>
	types::point<U> point() const
		{ return types::point<U>(_x, _y); }
};
			 
template <class T> class array_raster {
public:
	typedef T elem_type;
	typedef raster_iterator<array_raster<T>, T &> iterator;
	typedef const_raster_iterator<array_raster<T>, const T &> const_iterator;
private:
	typedef unsigned char byte;
	size_type	_width;
	size_type	_height;
	int			_pitch; //in bytes
	byte *		_data;
	std::vector<T *> _scans;
public:
	array_raster(void * data, size_type w, size_type h, order_type o = top_down) 
		: _data((byte*)data), _width(w), _height(h), _pitch(o == top_down ? sizeof(T) * (int)w : -sizeof(T) * (int)w)
		{ update_scans(); }
	// data is the first byte of row #0 (in the case of negative pitch it's NOT the first byte of array)
	array_raster(void * data, size_type w, size_type h, int pitch) 
		: _data((byte*)data), _width(w), _height(h), _pitch(pitch)
		{ update_scans(); }
	size_type width() const 
		{ return _width; }
	size_type height() const
		{ return _height; }
	size_extent dim() const
		{ return size_extent(_width, _height); }
	order_type order() const
		{ return _pitch > 0 : top_down : bottom_up; }
	const void * data() const
		{ return _data; }
	void * data() 
		{ return _data; }
	bool empty() const
		{ return _data == 0 || _width == 0 || _height == 0; }
	const T & operator () (size_type x, size_type y) const
		{ return _scans[y][x]; }
	const T & operator [] (point<size_type> p) const
		{ return _scans[p.y][p.x]; }
	const T & operator [] (point<int> p) const
		{ return _scans[p.y][p.x]; }
	T & operator () (size_type x, size_type y)
		{ return _scans[y][x]; }
	T & operator [] (point<size_type> p)
		{ return _scans[p.y][p.x]; }
	T & operator [] (point<int> p)
		{ return _scans[p.y][p.x]; }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, _height); }
	iterator begin()
		{ return iterator(* this, 0, 0); }
	iterator end()
		{ return iterator(* this, 0, _height); }
private:
	void update_scans() {
		_scans.resize(_height);
		unsigned int n = 0;
		typename std::vector<T *>::iterator it;
		for (it = _scans.begin(); it != _scans.end(); ++it)
			* it = reinterpret_cast<T *>(_data + (n++) * _pitch);
	}
};

template <class T> class const_array_raster {
public:
	typedef T elem_type;
	typedef const_raster_iterator<const_array_raster<T>, const T &> const_iterator;
private:
	size_type	_width;
	size_type	_height;
	order_type	_order;
	const T *	_data;
	std::vector<const T *> _scans;
public:
	const_array_raster(const T * data, size_type w, size_type h, order_type o = top_down) 
		: _data(data), _width(w), _height(h), _order(o)
	{
		_scans.resize(_height);
		unsigned int n = 0;
		if (_order == bottom_up) {
			typename std::vector<const T *>::reverse_iterator it;
			for (it = _scans.rbegin(); it != _scans.rend(); ++it)
				* it = reinterpret_cast<const T *>(_data + (n++) * _width);
		} else {
			typename std::vector<const T *>::iterator it;
			for (it = _scans.begin(); it != _scans.end(); ++it)
				* it = reinterpret_cast<const T *>(_data + (n++) * _width);
		}
	}
	size_type width() const 
		{ return _width; }
	size_type height() const
		{ return _height; }
	size_extent dim() const
		{ return size_extent(_width, _height); }
	order_type order() const
		{ return _order; }
	const T * data() const
		{ return _data; }
	bool empty() const
		{ return _data == 0 || _width == 0 || _height == 0; }
	const T & operator () (size_type x, size_type y) const
		{ return _scans[y][x]; }
	const T & operator [] (point<size_type> p) const
		{ return _scans[p.y][p.x]; }
	const T & operator [] (point<int> p) const
		{ return _scans[p.y][p.x]; }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, _height); }
};

template <class T>
array_raster<T> raster(T * data, size_type w, size_type h, order_type o = top_down) 
	{ return array_raster<T>(data, w, h, o); }
template <class T>
array_raster<T> raster(T * data, size_extent e, order_type o = top_down) 
	{ return array_raster<T>(data, e.sx, e.sy, o); }

template <class T>
const_array_raster<T> raster(const T * data, size_type w, size_type h, order_type o = top_down)
	{ return const_array_raster<T>(data, w, h, o); }
template <class T>
const_array_raster<T> raster(const T * data, size_extent e, order_type o = top_down)
	{ return const_array_raster<T>(data, e.sx, e.sy, o); }

//acts like const_raster, but always returns one and the same value from all positions
template <class T> class uniform_raster {
public:
	typedef T elem_type;
	typedef const_raster_iterator<uniform_raster<T>, const T &> const_iterator;
private:
	size_type	_width;
	size_type	_height;
	T			_value;
public:
	uniform_raster(const T &value, size_type w, size_type h) 
		: _value(value), _width(w), _height(h) {}
	uniform_raster(const T &value, size_extent e) 
		: _value(value), _width(e.sx), _height(e.sy) {}
	size_type width() const 
		{ return _width; }
	size_type height() const
		{ return _height; }
	size_extent dim() const
		{ return size_extent(_width, _height); }
	const T & value() const
		{ return _value; }
	bool empty() const
		{ return _width == 0 || _height == 0; }
	const T & operator () (size_type x, size_type y) const
		{ return _value; }
	const T & operator [] (point<size_type> p) const
		{ return _value; }
	const T & operator [] (point<int> p) const
		{ return _value; }
	T operator () (float x, float y) const
		{ return _value; }
	T operator [] (point<float> p) const
		{ return _value; }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, _height); }
};

template <class T>
inline uniform_raster<T> make_uniform_raster(const T & value, size_type w, size_type h) 
	{ return uniform_raster<T>(value, w, h); }

template <class T>
inline uniform_raster<T> make_uniform_raster(const T & value, size_extent e) 
	{ return uniform_raster<T>(value, e); }

template <class R, class OP, class T> class modify_raster {
public:
	typedef T elem_type;
	typedef const_raster_iterator<modify_raster<R,OP,T>, T> const_iterator;
private:
	const R & _raster;
	OP _op;
public:
	modify_raster(const R & raster, OP op)
		: _raster(raster), _op(op) { }
	size_type width() const
		{ return _raster.width(); }
	size_type height() const
		{ return _raster.height(); }
	size_extent dim() const
		{ return size_extent(_raster.width(), _raster.height()); }
	bool empty() const
		{ return _raster.empty(); }
	T operator () (size_type x, size_type y) const
		{ return _op(_raster(x, y)); }
	T operator [] (point<size_type> p) const
		{ return _op(_raster[p]); }
	T operator [] (point<int> p) const
		{ return _op(_raster[p]); }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, height()); }
};

template <class R, class OP>
modify_raster<R, OP &, typename OP::result_type> modify(const R & raster, OP & op)
	{ return modify_raster<R, OP &, typename OP::result_type>(raster, op); }

template <class R, class A, class T>
modify_raster<R, T (*)(A), T> modify(const R & raster, T (* fun)(A)) 
	{ return modify_raster<R, T (*)(A), T>(raster, fun); }

template <class R, class OP, class T> class transform_raster {
public:
	typedef T elem_type;
	typedef const_raster_iterator<transform_raster<R,OP,T>, T> const_iterator;
private:
	const R & _raster;
	OP _op;
public:
	transform_raster(const R & raster, OP op)
		: _raster(raster), _op(op) { }
	size_type size() const
		{ return _raster.size(); }
	size_type width() const
		{ return _raster.width(); }
	size_type height() const
		{ return _raster.height(); }
	size_extent dim() const
		{ return size_extent(_raster.width(), _raster.height()); }
	bool empty() const
		{ return _raster.empty(); }
	T operator () (size_type x, size_type y) const
		{ return _op(_raster, x, y); }
	T operator [] (point<size_type> p) const
		{ return _op(_raster, p.x, p.y); }
	T operator [] (point<int> p) const
		{ return _op(_raster, p.x, p.y); }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, height()); }
};

template <class R, class OP>
inline transform_raster<R, OP &, typename OP::result_type> transform(const R & raster, OP & op)
	{ return transform_raster<R, OP &, typename OP::result_type>(raster, op); }

template <class R, class I, class T>
inline transform_raster<R, T (*)(R, I, I), T> transform(const R & raster, T (* fun)(R, I, I)) 
	{ return transform_raster<R, T (*)(R, I, I), T>(raster, fun); }

template <class R1, class R2, class Op>
void transform(R1 & r1, const size_rect & where,
		   const R2 & r2, const size_point & from, Op & op)
{
	if (where.x1 >= r1.width() || where.y1 >= r1.height())
		return;
	if (from.x >= r2.width() || from.y >= r2.height())
		return;
	size_point bottom_right;
	bottom_right.x = std::min(
		std::min(where.x2, r1.width()), 
		where.x1 + (r2.width() - from.x));
	bottom_right.y = std::min(
		std::min(where.y2, r1.height()), 
		where.y1 + (r2.height() - from.y));
	size_point to = where.p1();
	size_point p;
	for (p.y = 0; p.y + where.y1 != bottom_right.y; ++p.y) {
		for (p.x = 0; p.x + where.x1 != bottom_right.x; ++p.x) {
			op(r1[to + p], r2[from + p]);
		}
	}
}

template <class R1, class R2, class Op>
inline void transform(R1 & r1, const R2 & r2, const size_point & from, Op & op)
	{ transform(r1, size_rect(0,0,r1.width(),r1.height()), r2, from, op); }

template <class R1, class R2, class Op>
inline void transform(R1 & r1, const R2 & r2, Op & op)
	{ transform(r1, r2, size_point(0,0), op); }

struct copy_op {
	template <class C1, class C2>
	void operator () (C1 & c1, const C2 & c2) const
		{ c1 = c2; }
};

template <class R1, class R2>
inline void copy(R1 & r1, const size_rect & where, const R2 & r2, const size_point & from = size_point(0,0))
	{ transform(r1, where, r2, from, copy_op()); }

template <class R1, class R2>
inline void copy(R1 & r1, const R2 & r2, const size_point & from)
	{ copy(r1, size_rect(0,0,r1.width(),r1.height()), r2, from); }

template <class R1, class R2>
inline void copy(R1 & r1, const R2 & r2)
	{ copy(r1, r2, size_point(0,0)); }

//vc71 workaround
template <class T, class E>
struct element_pointer_typedef
{
	typedef T (E::* pmem_type);
};

template <class R, class T> class projection_raster {
public:
	typedef T elem_type;
	typedef typename R::elem_type E;
	typedef typename element_pointer_typedef<T,E>::pmem_type pmem_type;
	typedef raster_iterator<projection_raster<R,T>, T &> iterator;
	typedef const_raster_iterator<projection_raster<R,T>, const T &> const_iterator;
private:
	R &			_raster;
	pmem_type	_pmem;
public:
	projection_raster(R & raster, pmem_type pmem)
		: _raster(raster), _pmem(pmem) { }
	size_type width() const
		{ return _raster.width(); }
	size_type height() const
		{ return _raster.height(); }
	size_extent dim() const
		{ return size_extent(_raster.width(), _raster.height()); }
	bool empty() const
		{ return _raster.empty(); }
	T & operator () (size_type x, size_type y)
		{ return _raster(x, y).*_pmem; }
	T & operator [] (point<size_type> p)
		{ return _raster[p].*_pmem; }
	T & operator [] (point<int> p)
		{ return _raster[p].*_pmem; }
	const T & operator () (size_type x, size_type y) const
		{ return _raster(x, y).*_pmem; }
	const T & operator [] (point<size_type> p) const
		{ return _raster[p].*_pmem; }
	const T & operator [] (point<int> p) const
		{ return _raster[p].*_pmem; }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, height()); }
	iterator begin()
		{ return iterator(* this, 0, 0); }
	iterator end()
		{ return iterator(* this, 0, height()); }
	template <class S> 
	projection_raster & operator = (const S & raster) {
		assert(width() == raster.width() && height() == raster.height());
		iterator i;
		for (i = begin(); i != end(); ++i)
			*i = raster[i];
		return * this;
	}
};

template <class R, class T, class E> //last argument to vc71 workaround
inline projection_raster<R, T> projection(R & raster, T (E::* pmem))
	{ return projection_raster<R, T>(raster, pmem); }

template <class R, class T> class const_projection_raster {
public:
	typedef T elem_type;
	typedef typename R::elem_type E;
	typedef typename element_pointer_typedef<T,E>::pmem_type pmem_type;
	typedef raster_iterator<const_projection_raster<R,T>, T &> iterator;
	typedef const_raster_iterator<const_projection_raster<R,T>, const T &> const_iterator;
private:
	const R &	_raster;
	pmem_type	_pmem;
public:
	const_projection_raster(const R & raster, pmem_type pmem)
		: _raster(raster), _pmem(pmem) { }
	size_type width() const
		{ return _raster.width(); }
	size_type height() const
		{ return _raster.height(); }
	size_extent dim() const
		{ return size_extent(_raster.width(), _raster.height()); }
	bool empty() const
		{ return _raster.empty(); }
	const T & operator () (size_type x, size_type y) const
		{ return _raster(x, y).*_pmem; }
	const T & operator [] (point<size_type> p) const
		{ return _raster[p].*_pmem; }
	const T & operator [] (point<int> p) const
		{ return _raster[p].*_pmem; }
	T operator () (float x, float y) const
		{ return interpolate(* this, x, y); }
	T operator [] (point<float> p) const
		{ return interpolate(* this, p.x, p.y); }
	const_iterator begin() const
		{ return const_iterator(* this, 0, 0); }
	const_iterator end() const
		{ return const_iterator(* this, 0, height()); }
};

template <class R, class T, class E> //last argument to vc71 workaround
inline const_projection_raster<R, T> const_projection(const R & raster, T (E::* pmem))
	{ return const_projection_raster<R, T>(raster, pmem); }

template <class R, class OP> void for_each(R & r, OP & op) {
	size_type x, y, w = r.width(), h = r.height();
	for (y = 0; y < h; ++y)
		for (x = 0; x < w; ++x)
			op(r, x, y);
}

template <class R>
class set_element {
public:
	typedef typename R::elem_type T;
	set_element(const T &v) : _v(v) {}
	void operator() (R &r, size_type x, size_type y)
		{ r(x,y) = _v; }
private:
	T _v;
};

template <class R> 
inline void fill(R & r, const typename R::elem_type &v) {
	for_each(r, set_element<R>(v) );
}

template <class R>
inline bool equal(const R & r1, const R & r2)
{
	if (r1.dim() != r2.dim())
		return false;
	return std::equal(r1.begin(), r1.end(), r2.begin());
}

} // img

#endif // IMGLIB_RASTER_H
